Shot Noise Analysis

Setup

Using R version 4.5.0 (2025-04-11) , tidyverse version 2.0.0, and isoorbi version 1.4.1.1.

# load packages
library(isoorbi) # for orbitrap function
library(dplyr) # for mutating data frames
library(forcats) # for recoding factors
library(ggplot2) # for data visualization
library(cowplot) # arrange multipanel plots

Trifluoroacetate (TFA)

Data from RAW (recommended)

# read files
raw_files <- 
  "data/shot_noise/tfa" |>
  orbi_find_raw() |>
  orbi_read_raw()

✔ [10.4s] orbi_read_raw() read TFA_1E4_AGC.raw (31.04 MB)

✔ [2.7s] orbi_read_raw() read TFA_M0_1uscan_120kRes.raw (6.62 MB)

✔ [11.4s] orbi_read_raw() read TFA_M0_1uscan_15kRes.raw (35.69 MB)

✔ [24.6s] orbi_read_raw() finished reading 3 files

# aggregate data
raw_agg <- raw_files |> orbi_aggregate_raw()

✔ [207ms] aggregate_files() aggregated file_info (3), scans (17.51 k), peaks
(91.67 k), and spectra (0) from 3 files

  → file_info: uid, file_path, RAW file, RAW file version, Creation date,
  Operator, Number of instruments, Description, Instrument model, Instrument
  name, Instrument method, Serial number, Software version, Firmware version,
  Units, Mass resolution, Number of scans, Number of ms2 scans, …, User text 3,
  and User text 4

  → scans: uid, scan, time.min, tic, it.ms, resolution, basePeakIntensity,
  rawOvFtT, intensCompFactor, agc, agcTarget, microscans,
  numberLockmassesFound, and analyzerTemperature

  → peaks: uid, scan, mzMeasured, intensity, and peakNoise

  → spectra: uid, scan, mz, and intensity

  → problems: uid, type, call, message, and condition

# fetch dataset of interest and calculate quants
ds <- 
  raw_agg |>
  orbi_get_data(
    file_info = c("uid", "filepath" = "file_path"),
    scans = c("uid", "scan", "time.min", "tic", "it.ms", "resolution", "microscans"),
    peaks = everything()
  ) |>
  rename("filename" = "uid", "scan.no" = "scan") |> 
  mutate(filename = gsub("\\.raw", "", filename)) |>
  cross_join(
    readr::read_tsv("data/shot_noise/tfa/isotopologs.tsv",
             show_col_types = FALSE) |>
      rename("compound" = "#Compound", "isotopocule" = "isotopolog")
  ) |>
  filter(abs(mass - mzMeasured) * 1000 < `Tolerance [mmu]`) |>
  # calculate ions incremental
  # todo: turn into a function (see Caj's notes on the papers)
  mutate(
    ions.incremental = intensity / peakNoise * 3.0 / z *
      sqrt(240000 / resolution) * sqrt(microscans),
    .after = "intensity"
  )

✔ [84ms] orbi_get_data() retrieved 91.67 k records from the combination of
file_info (3), scans (17.51 k), and peaks (91.67 k) via uid and scan

# data tfa
data_tfa <-
  ds |>
  orbi_flag_satellite_peaks() |>
  orbi_flag_weak_isotopocules(min_percent = 10) |>
  orbi_flag_outliers(agc_fold_cutoff = 2) |>
  orbi_define_basepeak("M0")

✔ [495ms] orbi_flag_satellite_peaks() flagged 5/52401 peaks in 1 isotopocules
("18O") as satellite peaks (0.0%)

ℹ column isotopocule was turned into a factor
✔ [13ms] orbi_flag_weak_isotopocules() confirmed there are no weak
isotopocules: all are detected in at least 10% of scans in each of the 9 data
groups (based on filename, compound, and isotopocule)

✔ [56ms] orbi_flag_outliers() confirmed that none of the 17513 scans are
outliers based on 2 fold AGC cutoff, i.e. based on scans below 1/2 and above 2
times the average number of ions tic * it.ms in the Orbitrap analyzer, in 3
data groups (based on filename)

✔ [285ms] orbi_define_basepeak() set M0 as the ratio denominator and calculated
34888 ratios for 2 isotopocules/base peak (13C and 18O)

shot_noise_tfa <- 
  data_tfa |>
  orbi_analyze_shot_noise() |>
  orbi_export_data_to_excel("output/shot_noise_tfa.xlsx")

✔ [20ms] orbi_analyze_shot_noise() analyzed the shot noise for 34883 peaks
(excluding 5 flagged peaks)
✔ [2.8s] orbi_export_data_to_excel() exported dataset (34883 rows, 27 columns)
to output/shot_noise_tfa.xlsx

Data from ISOX (obsolete)

# originally, this is how the data was read in
data_tfa <-
  "data/shot_noise/tfa" |>
  orbi_find_isox() |>
  orbi_read_isox() |>
  orbi_simplify_isox() |>
  orbi_flag_satellite_peaks() |>
  orbi_flag_weak_isotopocules(min_percent = 10) |>
  orbi_flag_outliers(agc_fold_cutoff = 2) |>
  orbi_define_basepeak("M0")

shot_noise_tfa <- 
  data_tfa |>
  orbi_analyze_shot_noise() |>
  orbi_export_data_to_excel("output/shot_noise_tfa.xlsx")

Bonus Figure: TFA satellite peaks

data_tfa |> orbi_plot_satellite_peaks()

Figure 1: Satellite peaks

Bonus Figure: TFA isotopocule coverage

data_tfa |> orbi_plot_isotopocule_coverage()

Figure 2: Isotopocule coverage

Bonus Figure: TFA ratios

data_tfa |> orbi_plot_raw_data(y = ratio)

Figure 3: Isotopocule ratios vs. M0

Figure 10: TFA shotnoise vs counts/time

# individual plots
tfa_vs_ions <-
  shot_noise_tfa |>
  orbi_filter_isox(filename = "TFA_M0_1uscan_15kRes") |>
  orbi_plot_shot_noise(x = "n_effective_ions") +
  labs(title = "vs counts") +
  theme(plot.title = element_text(hjust = 0.5))

tfa_vs_time <-
  shot_noise_tfa |>
  orbi_filter_isox(filename = "TFA_M0_1uscan_15kRes") |>
  orbi_plot_shot_noise(permil_target = 1) +
  labs(title = "vs time") +
  theme(plot.title = element_text(hjust = 0.5))

# combine
plot_grid(
  tfa_vs_ions + theme(legend.position.inside = c(0.85, 0.74)), 
  tfa_vs_time + theme(legend.position = "none"),
  align = "h", nrow = 1, axis = "tb"
)

Figure 4: Shot noise vs counts and vs time.

Extended Data Figure 6: TFA shotnoise at different IT

shot_noise_tfa |>
  orbi_filter_isox(filename = c("TFA_1E4_AGC", "TFA_M0_1uscan_15kRes")) |>
  # change color legend
  mutate(
    IT_info = filename |>
      fct_recode(
        "IT = 0.03 ms" = "TFA_1E4_AGC",
        "IT = 0.75 ms" = "TFA_M0_1uscan_15kRes"
      )
  ) |>
  orbi_plot_shot_noise(color = "IT_info") +
  # wrap by the ratio label
  facet_wrap(~ratio_label) +
  theme(legend.position = c(0.41, 0.75))

Figure 5: Shot noise at different ITs.

Extended Data Figure 7: TFA shotnoise at different resolutions

shot_noise_tfa |>
  orbi_filter_isox(filename = c("TFA_M0_1uscan_120kRes", "TFA_M0_1uscan_15kRes")) |>
  # change color legend
  mutate(
    res_info = filename |>
      fct_recode(
        "res. = 120k" = "TFA_M0_1uscan_120kRes",
        "res. = 15k" = "TFA_M0_1uscan_15kRes"
      )
  ) |>
  orbi_plot_shot_noise(color = "res_info") +
  # wrap by the ratio label
  facet_wrap(~ratio_label) +
  theme(legend.position = c(0.42, 0.75))

Figure 6: Shot noise at different resolutions.

Model Peptide (MRFA)

Data

# load, process, and export mrfa data
data_aas <- 
  "data/shot_noise/mrfa" |>
  orbi_find_isox() |>
  orbi_read_isox() |>
  orbi_simplify_isox() |>
  orbi_flag_satellite_peaks() |>
  orbi_flag_weak_isotopocules(min_percent = 10) |> 
  orbi_flag_outliers(agc_fold_cutoff = 2) |> 
  orbi_define_basepeak("M0") |>
  mutate(
    compound = compound |>
      fct_recode(
        "Alanine" = "Ala_imm",
        "Arginine" = "Arg_imm",
        "Methionine" = "Met_imm",
        "Phenylalanine" = "Phe_imm"
      )
  )

✔ [12ms] orbi_read_isox() loaded 3853 peaks for 4 compounds (Ala_imm, Arg_imm,
Met_imm, and Phe_imm) with 6 isotopocules (M0, 15N, 13C, 2H, 33S, and 34S) from
MRFA_M1_40NCE_10min_10uscans.isox

✔ [40ms] orbi_read_isox() read .isox data from 1 file

✔ [2ms] orbi_simplify_isox() kept columns filepath, filename, scan.no,
time.min, compound, isotopocule, ions.incremental, tic, and it.ms

✔ [44ms] orbi_flag_satellite_peaks() flagged 59/3853 peaks in 3 isotopocules
("M0", "15N", "13C") as satellite peaks (1.5%)

✔ [19ms] orbi_flag_weak_isotopocules() confirmed there are no weak
isotopocules: all are detected in at least 10% of scans in each of the 18 data
groups (based on filename, compound, and isotopocule)

✔ [9ms] orbi_flag_outliers() confirmed that none of the 230 scans are outliers
based on 2 fold AGC cutoff, i.e. based on scans below 1/2 and above 2 times the
average number of ions tic * it.ms in the Orbitrap analyzer

✔ [21ms] orbi_define_basepeak() set M0 as the ratio denominator and calculated
2915 ratios for 5 isotopocules/base peak (15N, 13C, 2H, 33S, and 34S)

shot_noise_aas <-
  data_aas |>
  orbi_analyze_shot_noise() |> 
  orbi_export_data_to_excel("output/shot_noise_MRFA.xlsx") 

✔ [16ms] orbi_analyze_shot_noise() analyzed the shot noise for 2874 peaks
(excluding 41 flagged peaks)

✔ [227ms] orbi_export_data_to_excel() exported dataset (2874 rows, 19 columns)
to output/shot_noise_MRFA.xlsx

Table

# example of first few rows
shot_noise_aas |>
  arrange(compound, isotopocule, scan.no) |>
  select(compound, scan.no, time.min, isotopocule,
         ratio, ratio_rel_se.permil, shot_noise.permil) |>
  head(10) |>
  knitr::kable()

compound	scan.no	time.min	isotopocule	ratio	ratio_rel_se.permil	shot_noise.permil
Alanine	1	0.059	15N	0.0062453	NaN	73.71127
Alanine	2	0.102	15N	0.0060919	12.43207	52.42764
Alanine	3	0.146	15N	0.0068807	37.74616	42.08467
Alanine	4	0.189	15N	0.0059526	32.61878	36.78579
Alanine	5	0.233	15N	0.0060549	26.57095	32.90351
Alanine	6	0.277	15N	0.0055040	29.99353	30.32862
Alanine	7	0.320	15N	0.0060244	25.50106	28.18264
Alanine	8	0.364	15N	0.0065182	23.43771	26.31067
Alanine	9	0.408	15N	0.0064846	21.36449	24.76418
Alanine	10	0.451	15N	0.0061980	19.10467	23.50447

Bonus Figure: Amino Acids satellite peaks

data_aas |> orbi_plot_satellite_peaks()

Figure 7: Satellite peaks

Bonus Figure: Amino Acids isotopocule coverage

data_aas |> orbi_plot_isotopocule_coverage()

Figure 8: Isotopocule coverage

Bonus Figure: Amino Acids ratios

data_aas |> orbi_plot_raw_data(y = ratio)

Figure 9: Isotopocule ratios vs. M0

Extended Data Figure 8: Amino Acids Shot Noise

shot_noise_aas |>
  filter(compound == "Methionine" | isotopocule != "2H") |>
  orbi_plot_shot_noise()

Figure 10: amino acids shotnoise